library(MVN)
library(haven)

# Load Data ---------------------------------------------------------------
setwd("./Code/")
s2 <- read_dta("../data/study2_raw.dta")
source("./main/study3_cleaning.R")
source("./main/study4_cleaning.R")


# Study 2 -----------------------------------------------------------------
# Merging grid presentation and recoding to pro-trait direction
s2$populism_for_1[which(is.na(s2$populism_for_1))] <- s2$populism_full_1[which(is.na(s2$populism_for_1))]
s2$populism_for_2[which(is.na(s2$populism_for_2))] <- s2$populism_full_2[which(is.na(s2$populism_for_2))]
s2$populism_for_3[which(is.na(s2$populism_for_3))] <- s2$populism_full_3[which(is.na(s2$populism_for_3))]
s2$populism_rev_1[which(is.na(s2$populism_rev_1))] <- s2$populism_full_4[which(is.na(s2$populism_rev_1))]
s2$populism_rev_2[which(is.na(s2$populism_rev_2))] <- s2$populism_full_5[which(is.na(s2$populism_rev_2))]
s2$populism_rev_3[which(is.na(s2$populism_rev_3))] <- s2$populism_full_6[which(is.na(s2$populism_rev_3))]
s2$populism_rev_1 <- s2$populism_rev_1*-1 + 6
s2$populism_rev_2 <- s2$populism_rev_2*-1 + 6
s2$populism_rev_3 <- s2$populism_rev_3*-1 + 6


s2$hs_forward_1[which(is.na(s2$hs_forward_1))] <- s2$hs_full_1[which(is.na(s2$hs_forward_1))]
s2$hs_forward_2[which(is.na(s2$hs_forward_2))] <- s2$hs_full_2[which(is.na(s2$hs_forward_2))]
s2$hs_forward_3[which(is.na(s2$hs_forward_3))] <- s2$hs_full_3[which(is.na(s2$hs_forward_3))]
s2$hs_reversed_1[which(is.na(s2$hs_reversed_1))] <- s2$hs_full_4[which(is.na(s2$hs_reversed_1))]
s2$hs_reversed_2[which(is.na(s2$hs_reversed_2))] <- s2$hs_full_5[which(is.na(s2$hs_reversed_2))]
s2$hs_reversed_3[which(is.na(s2$hs_reversed_3))] <- s2$hs_full_6[which(is.na(s2$hs_reversed_3))]
s2$hs_reversed_1 <- s2$hs_reversed_1*-1 + 6
s2$hs_reversed_2 <- s2$hs_reversed_2*-1 + 6
s2$hs_reversed_3 <- s2$hs_reversed_3*-1 + 6

v_antidem <- grep("antidem", names(s2), ignore.case = T, value = T, fixed = F)
v_conspre <- grep("conspre", names(s2), ignore.case = T, value = T, fixed = F)
v_pop <- grep("populism_", names(s2), ignore.case = T, value = T, fixed = F)
v_pop <- grep("full", v_pop, ignore.case = T, value = T, fixed = F, invert = T)
v_hs <- grep("hs_", names(s2), ignore.case = T, value = T, fixed = F)
v_hs <- grep("full", v_hs, ignore.case = T, value = T, fixed = F, invert = T)


s2_desc <- data.frame(scale = c("Anti-Dem", "Conspiracy", "Populism", "HS"),
                      skew_med = c(abs(median(mvn(data = s2[, v_antidem], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = s2[, v_conspre], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = s2[, v_pop], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = s2[, v_hs], mvn_test = "mardia")$descriptives[, "Skew"]))),
                      skew_range = c(paste(abs(range(mvn(data = s2[, v_antidem], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = s2[, v_conspre], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = s2[, v_pop], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = s2[, v_hs], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-")),
                      kurt_med = c(abs(median(mvn(data = s2[, v_antidem], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = s2[, v_conspre], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = s2[, v_pop], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = s2[, v_hs], mvn_test = "mardia")$descriptives[, "Kurtosis"]))),
                      kurt_range = c(paste(abs(range(mvn(data = s2[, v_antidem], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = s2[, v_conspre], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = s2[, v_pop], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = s2[, v_hs], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"))
                      )

# write.csv(s2_desc, file = "s2_mvn.csv", row.names = F)


mvn(data = s2[, v_antidem], mvn_test = "mardia")[1:3]
mvn(data = s2[, v_conspre], mvn_test = "mardia")[1:3]
mvn(data = s2[, v_pop], mvn_test = "mardia")[1:3]
mvn(data = s2[, v_hs], mvn_test = "mardia")[1:3]

# Study 3 -----------------------------------------------------------------
v <- grep("nfc_", names(mt), value = T)
v_nfc <- v[1:15]

v <- grep("pop_", names(mt), value = T)
v_pop <- v[1:6]

v <- grep("viol_", names(mt), value = T)
v_viol <- v[1:8]

v <- grep("consp_", names(mt), value = T)
v_consp <- v[1:9]


v <- grep("rr_", names(mt), value = T)
v_rr <- v[1:4]

v <- grep("hs_", names(mt), value = T)
v_hs <- v[1:6]


s3_desc <- data.frame(scale = c("NFC", "Conspiracy", "Populism", "HS", "Violence", "RR"),
                      skew_med = c(abs(median(mvn(data = mt[, v_nfc], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = mt[, v_consp], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = mt[, v_pop], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = mt[, v_hs], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = mt[, v_viol], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = mt[, v_rr], mvn_test = "mardia")$descriptives[, "Skew"]))),
                      skew_range = c(paste(abs(range(mvn(data = mt[, v_nfc], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_consp], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_pop], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_hs], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_viol], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_rr], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-")
                                  ),
                      kurt_med = c(abs(median(mvn(data = mt[, v_nfc], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = mt[, v_consp], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = mt[, v_pop], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = mt[, v_hs], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = mt[, v_viol], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = mt[, v_rr], mvn_test = "mardia")$descriptives[, "Kurtosis"]))),
                      kurt_range = c(paste(abs(range(mvn(data = mt[, v_nfc], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_consp], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_pop], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_hs], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_viol], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = mt[, v_rr], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-")
                                     )
)

write.csv(s3_desc, file = "s3_mvn.csv", row.names = F)

mvn(data = mt[, v_nfc], mvn_test = "mardia")[1:3]
mvn(data = mt[, v_pop], mvn_test = "mardia")[1:3]
mvn(data = mt[, v_viol], mvn_test = "mardia")[1:3]
mvn(data = mt[, v_consp], mvn_test = "mardia")[1:3]
mvn(data = mt[, v_rr], mvn_test = "mardia")[1:3]
mvn(data = mt[, v_hs], mvn_test = "mardia")[1:3]

# Study 4 -----------------------------------------------------------------
v <- grep("nfc_", names(bov), value = T)
v_nfc <- v[1:15]

v <- grep("pop_", names(bov), value = T)
v_pop <- v[2:7]

v <- grep("viol_", names(bov), value = T)
v_viol <- v[1:8]

v <- grep("consp_", names(bov), value = T)
v_consp <- v[1:9]

v <- grep("rr_", names(bov), value = T)
v_rr <- v[1:4]

v <- grep("hs_", names(bov), value = T)
v_hs <- v[1:6]

v <- grep("antidemoc_", names(bov), value = T)
v_antidem <- v[1:8]

s4_desc <- data.frame(scale = c("NFC", "Conspiracy", "Populism", "HS", "Violence", "RR", "Anti-Dem"),
                      skew_med = c(abs(median(mvn(data = bov[, v_nfc], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = bov[, v_consp], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = bov[, v_pop], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = bov[, v_hs], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = bov[, v_viol], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = bov[, v_rr], mvn_test = "mardia")$descriptives[, "Skew"])),
                                   abs(median(mvn(data = bov[, v_antidem], mvn_test = "mardia")$descriptives[, "Skew"]))),
                      skew_range = c(paste(abs(range(mvn(data = bov[, v_nfc], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_consp], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_pop], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_hs], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_viol], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_rr], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_antidem], mvn_test = "mardia")$descriptives[, "Skew"])), collapse = "-")),
                      kurt_med = c(abs(median(mvn(data = bov[, v_nfc], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = bov[, v_consp], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = bov[, v_pop], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = bov[, v_hs], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = bov[, v_viol], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = bov[, v_rr], mvn_test = "mardia")$descriptives[, "Kurtosis"])),
                                   abs(median(mvn(data = bov[, v_antidem], mvn_test = "mardia")$descriptives[, "Kurtosis"]))),
                      kurt_range = c(paste(abs(range(mvn(data = bov[, v_nfc], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_consp], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_pop], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_hs], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_viol], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_rr], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-"),
                                     paste(abs(range(mvn(data = bov[, v_antidem], mvn_test = "mardia")$descriptives[, "Kurtosis"])), collapse = "-")
                      )
)

write.csv(s4_desc, file = "s4_mvn.csv", row.names = F)

mvn(data = bov[, v_nfc], mvn_test = "mardia")[1:3]
mvn(data = bov[, v_pop], mvn_test = "mardia")[1:3]
mvn(data = bov[, v_viol], mvn_test = "mardia")[1:3]
mvn(data = bov[, v_consp], mvn_test = "mardia")[1:3]
mvn(data = bov[, v_rr], mvn_test = "mardia")[1:3]
mvn(data = bov[, v_hs], mvn_test = "mardia")[1:3]
mvn(data = bov[, v_antidem], mvn_test = "mardia")[1:3]
